//
// File: xgeqp3.cpp
//
// MATLAB Coder version            : 5.4
// C/C++ source code generated on  : 25-Apr-2025 12:42:26
//

// Include Files
#include "xgeqp3.h"
#include "freqDomainHRV_rtwutil.h"
#include "rt_nonfinite.h"
#include "xnrm2.h"
#include "coder_array.h"
#include <cmath>
#include <cstring>
#include <string.h>

// Function Definitions
//
// Arguments    : ::coder::array<double, 2U> &A
//                double tau_data[]
//                int *tau_size
//                int jpvt[2]
// Return Type  : void
//
namespace coder {
namespace internal {
namespace lapack {
void xgeqp3(::coder::array<double, 2U> &A, double tau_data[], int *tau_size,
            int jpvt[2])
{
  int m;
  boolean_T guard1{false};
  m = A.size(0);
  *tau_size = A.size(0);
  if (*tau_size > 2) {
    *tau_size = 2;
  }
  if (*tau_size - 1 >= 0) {
    std::memset(&tau_data[0], 0, *tau_size * sizeof(double));
  }
  guard1 = false;
  if (A.size(0) == 0) {
    guard1 = true;
  } else {
    int ix;
    ix = A.size(0);
    if (ix > 2) {
      ix = 2;
    }
    if (ix < 1) {
      guard1 = true;
    } else {
      double vn1[2];
      double vn2[2];
      double work[2];
      double temp;
      int ma;
      int minmn;
      ma = A.size(0);
      minmn = A.size(0);
      if (minmn > 2) {
        minmn = 2;
      }
      jpvt[0] = 1;
      work[0] = 0.0;
      temp = blas::xnrm2(m, A, 1);
      vn1[0] = temp;
      vn2[0] = temp;
      jpvt[1] = 2;
      work[1] = 0.0;
      temp = blas::xnrm2(m, A, ma + 1);
      vn1[1] = temp;
      vn2[1] = temp;
      for (int i{0}; i < minmn; i++) {
        double atmp;
        double beta1;
        int b_i;
        int ii;
        int ii_tmp;
        int ip1;
        int knt;
        int lastv;
        int mmi;
        int pvt;
        ip1 = i + 2;
        ii_tmp = i * ma;
        ii = ii_tmp + i;
        mmi = m - i;
        ix = 0;
        if ((2 - i > 1) && (std::abs(vn1[1]) > std::abs(vn1[i]))) {
          ix = 1;
        }
        pvt = i + ix;
        if (pvt != i) {
          ix = pvt * ma;
          for (lastv = 0; lastv < m; lastv++) {
            knt = ix + lastv;
            temp = A[knt];
            b_i = ii_tmp + lastv;
            A[knt] = A[b_i];
            A[b_i] = temp;
          }
          ix = jpvt[pvt];
          jpvt[pvt] = jpvt[i];
          jpvt[i] = ix;
          vn1[pvt] = vn1[i];
          vn2[pvt] = vn2[i];
        }
        if (i + 1 < m) {
          atmp = A[ii];
          ix = ii + 2;
          tau_data[i] = 0.0;
          if (mmi > 0) {
            temp = blas::xnrm2(mmi - 1, A, ii + 2);
            if (temp != 0.0) {
              beta1 = rt_hypotd_snf(A[ii], temp);
              if (A[ii] >= 0.0) {
                beta1 = -beta1;
              }
              if (std::abs(beta1) < 1.0020841800044864E-292) {
                knt = 0;
                b_i = ii + mmi;
                do {
                  knt++;
                  for (lastv = ix; lastv <= b_i; lastv++) {
                    A[lastv - 1] = 9.9792015476736E+291 * A[lastv - 1];
                  }
                  beta1 *= 9.9792015476736E+291;
                  atmp *= 9.9792015476736E+291;
                } while ((std::abs(beta1) < 1.0020841800044864E-292) &&
                         (knt < 20));
                beta1 = rt_hypotd_snf(atmp, blas::xnrm2(mmi - 1, A, ii + 2));
                if (atmp >= 0.0) {
                  beta1 = -beta1;
                }
                tau_data[i] = (beta1 - atmp) / beta1;
                temp = 1.0 / (atmp - beta1);
                for (lastv = ix; lastv <= b_i; lastv++) {
                  A[lastv - 1] = temp * A[lastv - 1];
                }
                for (lastv = 0; lastv < knt; lastv++) {
                  beta1 *= 1.0020841800044864E-292;
                }
                atmp = beta1;
              } else {
                tau_data[i] = (beta1 - A[ii]) / beta1;
                temp = 1.0 / (A[ii] - beta1);
                b_i = ii + mmi;
                for (lastv = ix; lastv <= b_i; lastv++) {
                  A[lastv - 1] = temp * A[lastv - 1];
                }
                atmp = beta1;
              }
            }
          }
          A[ii] = atmp;
        } else {
          tau_data[i] = 0.0;
        }
        if (i + 1 < 2) {
          int jA;
          atmp = A[ii];
          A[ii] = 1.0;
          jA = (ii + ma) + 1;
          if (tau_data[0] != 0.0) {
            lastv = mmi - 1;
            ix = (ii + mmi) - 1;
            while ((lastv + 1 > 0) && (A[ix] == 0.0)) {
              lastv--;
              ix--;
            }
            pvt = 1;
            ix = jA;
            int exitg1;
            do {
              exitg1 = 0;
              if (ix <= jA + lastv) {
                if (A[ix - 1] != 0.0) {
                  exitg1 = 1;
                } else {
                  ix++;
                }
              } else {
                pvt = 0;
                exitg1 = 1;
              }
            } while (exitg1 == 0);
          } else {
            lastv = -1;
            pvt = 0;
          }
          if (lastv + 1 > 0) {
            if (pvt != 0) {
              work[0] = 0.0;
              knt = 0;
              for (ii_tmp = jA; ma < 0 ? ii_tmp >= jA : ii_tmp <= jA;
                   ii_tmp += ma) {
                temp = 0.0;
                b_i = ii_tmp + lastv;
                for (ix = ii_tmp; ix <= b_i; ix++) {
                  temp += A[ix - 1] * A[(ii + ix) - ii_tmp];
                }
                work[knt] += temp;
                knt++;
              }
            }
            if (!(-tau_data[0] == 0.0)) {
              for (ii_tmp = 0; ii_tmp < pvt; ii_tmp++) {
                if (work[0] != 0.0) {
                  temp = work[0] * -tau_data[0];
                  b_i = lastv + jA;
                  for (knt = jA; knt <= b_i; knt++) {
                    A[knt - 1] = A[knt - 1] + A[(ii + knt) - jA] * temp;
                  }
                }
                jA += ma;
              }
            }
          }
          A[ii] = atmp;
        }
        for (ii_tmp = ip1; ii_tmp < 3; ii_tmp++) {
          ix = i + ma;
          if (vn1[1] != 0.0) {
            temp = std::abs(A[ix]) / vn1[1];
            temp = 1.0 - temp * temp;
            if (temp < 0.0) {
              temp = 0.0;
            }
            beta1 = vn1[1] / vn2[1];
            beta1 = temp * (beta1 * beta1);
            if (beta1 <= 1.4901161193847656E-8) {
              if (i + 1 < m) {
                temp = blas::xnrm2(mmi - 1, A, ix + 2);
                vn1[1] = temp;
                vn2[1] = temp;
              } else {
                vn1[1] = 0.0;
                vn2[1] = 0.0;
              }
            } else {
              vn1[1] *= std::sqrt(temp);
            }
          }
        }
      }
    }
  }
  if (guard1) {
    jpvt[0] = 1;
    jpvt[1] = 2;
  }
}

} // namespace lapack
} // namespace internal
} // namespace coder

//
// File trailer for xgeqp3.cpp
//
// [EOF]
//
